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Abstract 

We use renormalization group methods to derive equations of motion 
for large scale variables in fluid dynamics. The large scale variables are 
averages of the underlying continuum variables over cubic volumes, and 
naturally live on a lattice. The resulting lattice dynamics represents a 
perfect discretization of continuum physics, i.e. grid artifacts are com- 
pletely eliminated. Perfect equations of motion are derived for static, 
slow flows of incompressible, viscous fluids. For Hagen-Poiseuille flow 
in a channel with square cross section the equations reduce to a perfect 
discretization of the Poisson equation for the velocity fleld with Dirichlet 
boundary conditions. The perfect large scale Poisson equation is used 
in a numerical simulation, and is shown to represent the continuum flow 
exactly. For non-square cross sections we use a numerical iterative pro- 
cedure to derive flow equations that are approximately perfect. 
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1 Introduction 



Problems of fluid flow are of great practical importance in such areas as mechanical 
engineering and meteorology. A full theoretical understanding of fluid flow is very 
difficult because of the complex dynamics. In particular, in the regime of fully 
developed turbulence we face a dynamical system with fluctuations on all length 
scales. The small scale structures influence the dynamics at larger scales in a way 
that is hard to control analytically. Numerical simulations are therefore an essential 
tool in fluid dynamics. However, with present day computers a full simulation 
of realistic problems from first principles is beyond reach. Already simple model 
problems exhaust the biggest machines available. Hence, in practical applications 
one makes additional assumptions, especially about the small scale dynamics, which 
are difficult to justify theoretically, and which involve free parameters. In fact, the 
primary systematic error in numerical simulations of fluid flow is due to the finite 
grid size that is necessary to discretize the continuum Navier-Stokes equations. 

The main purpose of this paper is to point out that perfect discretizations, 
i.e. ones that are completely free of grid artifacts, exist and can be constructed 
explicitly. The existence of perfect discretizations is a direct consequence of Wilson's 
renormalization group theory |I[ . The idea of the renormalization group is to deal 
with the fluctuations at different length scales step by step in scale. In this way 
the influence of small scale fluctuations on the large scale dynamics can, at least 
in principle, be controlled exactly. The renormalization group approach has been 
useful in the study turbulence and has been applied to various problems previously 
|T^, |]1T| for some examples). 



see 



Recently, however, in the context of lattice quantum field theories so-called per- 
fect actions have been constructed for various models using the renormalization 
group. Originally, Hasenfratz and Niedermayer iterated a renormalization group 
transformation to obtain a perfect discretization for the 2-dimensional 0(3) non- 
linear (T-model, and they demonstrated that grid artifacts were eliminated even on 
coarse lattices 0. The same method has been applied to lattice fermions 0, ^, to 
4-dimensional pure SU{3) gauge theory p and to full QCD [||. Here we use similar 
methods in classical physics. In principle, the renormalization group can be used 
to derive perfect discretizations for any differential equation. Examples of great 
practical importance are the Laplace and Poisson equations. Maxwell's equations 
and the Navier-Stokes equation. A perfect discretization of the Laplace equation 
was obtained by Bell and Wilson in an early study of scalar lattice field theory [^. 
The perfect action for a classical gauge field — and hence a perfect discretization 
of Maxwell's equations — has been constructed in refs.|^, ^. In this work we con- 
centrate on the Navier-Stokes equations. Here we do not treat the problem in its 
full complexity yet. Specifically, we restrict ourselves to incompressible fluids and 
we assume static, slow flows. 
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The first step in our program is to define large scale lattice variables from the 
underlying continuum fields. In the case of the Navier-Stokes equations these are 
the velocity and pressure fields. We average the continuum pressure over a cube 
of side length a to define the large scale pressure variables, which then naturally 
live on a lattice of spacing a consisting of the cube centers. Similarly, the velocity 
field is averaged over the interface between two neighboring cubes, and hence the 
corresponding lattice field lives on the links connecting neighboring lattice points. 
This construction ensures that the continuity equation maintains a simple form 
on the lattice. By construction we know exactly how the large scale variables are 
related to the continuum fields. The goal of the next step is to derive equations 
of motion for the large scale fields. Solving these equations yields exact results for 
the averaged continuum quantities. This is in contrast to standard discretization 
procedures, where the lattice field represents an approximation of the continuum 
field at the same point. The approximation becomes exact only in the continuum 
limit a — *■ 0, because the value of the field at a specific point is influenced by the 
dynamics at arbitrarily small scales. A perfect discretization, on the other hand, 
ensures exact results for averaged continuum quantities already at finite a. Still, the 
values of the continuum fields can be reconstructed from the lattice data by using 
so-called perfect fields. 

The second step is to derive exact equations of motion for the lattice fields. One 
possibility to do that is to iterate a discrete renormalization group transformation. 
This is done by starting on a very fine lattice (in fact, at the end the lattice spacing 
of the fine lattice is sent to zero), and then performing an infinite sequence of block- 
ing steps, in which the lattice fields are averaged over blocks, and effective equations 
of motion for the block variables are derived. Here we proceed differently by block- 
ing directly from the continuum. The result is identical to the discrete, iterative 
procedure, but the derivation is more transparent. Wilson was the first to discover 
this technique 0, which has not been very well known until recently. For example, 
blocking from the continuum was also used in the derivation of the QCD perfect 
action . In contrast to quantum field theory there is no action for the full Navier- 
Stokes problem. Still, blocking from the continuum can be directly applied to the 
equations of motion. The constraint that identifies the large scale lattice variables 
with the averaged continuum fields is implemented through a Lagrange multiplier. 
Solving the equations for the continuum variables in terms of the lattice fields and 
plugging the result back into the constraint equations leads to the equations for the 
large scale variables. 

For numerical applications it is essential that the lattice equations of motion are 
as local as possible. In this paper we distinguish three types of locality: ultralocal, 
local and nonlocal. A lattice equation is called ultralocal if it couples the lattice field 
at a given site to a finite number of neighboring sites only. Consider, for example. 
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the 1- dimensional Laplace equation 



0. 



(1.1) 



Its standard lattice discretization 



-[/(x + a)-2/(x) + /(x-a)] =0 



(1.2) 



is ultralocal because it couples the field at the point x to its nearest neighbors at 
x + a and x — a only. In general, a perfect discretization will not lead to an ultralo- 
cal coupling. The couplings to fields at distant sites are non-zero, but exponentially 
suppressed. In that case we say that the lattice equation is local. By optimizing 
the parameters of the renormalization group transformation one can control the 
strength of the exponential decay. It turns out that the parameters can be chosen 
such that the perfect equations become ultralocal in one dimension. In practice, even 
in higher dimensions this leads to extremely local perfect lattice equations, whose 
exponentially suppressed couplings to distant neighbors can safely be neglected. In 
particular, as it will become clear later, even the ultralocal eq.( [1.2|) turns out to 
be perfect. This would not be the case for the common approach to discretization 
problems, where one interprets lattice derivatives as approximations of continuum 
derivatives at the same point. Still, even then one can find a perfect discretiza- 
tion. It is a well known result in engineering that a smooth function f{x) can be 
reconstructed from its values f{na) on a lattice by 



Here we have assumed that f{x) is infinitely differentiable and that its Fourier 
transform exists. The above equation yields 



This exactly represents the continuum second derivative of the function at the lattice 
point ma. However, now the lattice derivative comprises terms at arbitrarily large 
distances, and their contributions are suppressed only power-like. This is what we 
call a nonlocal coupling. In this case the contributions from large distances cannot 
be neglected and therefore this equation is not of practical use. In contrast, the 
advantage of using block averaged fields is that their perfect discretizations are 
local. 

As opposed to previous applications of perfect discretizations, which dealt with 
infinite volumes or finite systems with periodic boundary conditions, the incorpo- 
ration of more general boundary conditions is essential in fiuid dynamics. While it 
is straightforward to implement periodic boundary conditions, it is more difficult to 
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treat fixed boundary conditions, which occur in typical fiuid dynamics applications. 
Still, it turns out that fixed boundary conditions can be represented exactly in the 
perfect lattice equations. Here we restrict ourselves to channels whose cross section 
is a rectangular polygon. Arbitrarily curved boundaries can in principle be handled 
using curvelinear coordinates, but we do not elaborate on this issue here. 

The paper is organized as follows. In section 2 we introduce the large scale lattice 
pressure and velocity variables by blocking from the continuum, and we investigate 
their behavior under discrete renormalization group transformations. Section 3 con- 
tains the derivation of perfect equations of motion for Hagen-Poiseuille fiow. The 
parameter in the renormalization group transformation is optimized for ultralocal- 
ity in two dimensions. We verify explicitly that also in three dimensions continuum 
physics is exactly reproduced. In section 4 we show analytically that the equations 
of motion indeed satisfy a factor 2 renormalization group step. Section 5 presents 
results of a numerical calculation of the lattice fields. In section 6 we construct 
the perfect velocity field that allows reconstruction of the continuum field from the 
lattice data. Section 7 deals with static, slow flows in general. We derive the perfect 
lattice version of the corresponding Navicr-Stokes equations. In section 8 we present 
a perturbative method of introducing the non-linear term to equation of motion. Fi- 
nally, section 9 contains concluding remarks on how to extend our construction to 
the full Navier-Stokes equations and to arbitrarily shaped boundaries. We also 
speculate about applications of perfect discretizations to numerical simulations of 
turbulent flows. 



2 Large scale lattice variables 

The continuum Navier-Stokes equations are formulated in terms of pressure and 
velocity variables. Here we consider static flows in d dimensions, and hence pressure 
p{y) and velocity Vi{y) with i e {1,2, ...,d} are functions of the position vector y 
only. For incompressible fluids the velocity fleld obeys the continuity equation 



To deflne the large scale lattice pressure variables the continuum pressure fleld 
p{y) is averaged over a d-dimensional cube Cx of side length a centered at x, 



The variables Px naturally live on a lattice of spacing a formed by the cube centers 
X. Similarly, for the velocity fleld we deflne 



diVi{y) = 0. 



(2.1) 




(2.2) 




(2.3) 
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where fi^x is the {d — l)-dimensional face separating the cubes c^_ai/2 ^^'^ '^x+ai/2- 
Note that i is the unit vector in the i-direction. It is natural to associate the velocity 
variables with the links connecting neighboring lattice sites x — ai/2 and x + ai/2. 
Now the z-component of x is a half-integer multiple of the lattice spacing. The 
lattice velocity variables can be interpreted as volume flow rates per cross sectional 
area of the faces separating neighboring cubes. This definition of the lattice velocity 
field is natural because the lattice continuity equation then assumes an ultralocal 
form. It is simply given by 



SVx 



hai/2 



^d-l E [J 



+ ai/2 



i,x — ai/2 



dc. 



d^y d,Vi{y) = 0. 



(2.4) 



Here we have used Gauss' law together with the continuum continuity equation. 
The definition of lattice variables by blocking from the continuum can be viewed as 
a renormalization group transformation with an infinite blocking factor. In fact, the 
same result could be obtained by iterating block factor 2 transformations starting 
from an arbitrarily fine lattice. To identify the corresponding factor 2 transformation 
we now block from the continuum to a lattice of twice the lattice spacing a' = 2a. 
The continuum pressure variables are then integrated over cubes c^/ of side length 
a' which contain 2^^ original cubes c^.. Hence, the resulting pressure variable 

i /, d'y piy) = hT.^,[d'y piy) = ^ E (2.5) 



p' 



2<^ — , a" 



on the coarse lattice is a block average of the pressure variables on the fine lattice. 
We use x G x' to denote that the cube belongs to the block c'^,. The corresponding 
geometry is illustrated in figure H. Similarly, for the velocity field we have 



P'(x') 



P(x) 



P(x+a i) 



V(x) 



V'(x') 



V(x-aj) 



P'(x'+a'i) 



Figure 1: The geometry of the block factor 2 renormalization group transformation. 



V , 



d'-'y Viiy) 
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= ^ E ^ /, d'~'y <y) = ^ E y^.■ (2-6) 

Here x & x' denotes that the sum includes faces fi^x on the fine lattice that belong to 
the face Z/ ^., on the coarse lattice. Of course, by construction the coarse lattice veloc- 
ity variables also obey the continuity equation. This block factor 2 renormalization 
group transformation (RGT) will be used later to demonstrate that the structure of 
the perfect equations of motion reproduces itself under renormalization. 



3 Perfect equations of motion for Hagen-Poiseuille 
flow 

For simplicity wc first discuss the derivation of perfect equations of motion in the 
context of Hagen-Poiseuille flow. The general equations for static, slow flows are 
derived in section 6. We start from the Navier-Stokes equation 

dtVi{y) + {vj{y)dj)vi{y) = ~dip{y) + vdjdjVi{y). (3.1) 

Here V is the kinematic viscosity of the fluid, and p is its density. We restrict 
ourselves to incompressible flows, i.e. to a constant density, and thus the continuity 
equation applies. The appropriate boundary condition for a viscous fluid at a flxed 
wall is Vi{y) — for all %. Furthermore, we consider static Hagen-Poiseuille flow, i.e. 
a constant pressure gradient 5p, in a channel with square cross section A — L^~^ 
along the rf-direction. Hence, the pressure is given by -piy) = Sp yd- The velocity 
is along the channel, i.e. only Vd{y) is nonzero and in addition independent of yd- 
Then the equations reduce to the {d — l)-dimensional Poisson equation 

didMy) = — , (3.2) 
up 

with Dirichlet boundary conditions Vd{y) = at the edge of the square. This 
equation follows from a variational principle of an action 

= / d'-'y[ldMy)dMy) + -Vd{y)]. (3.3) 

Ja 2 up 

Note that here i runs from 1 to (o? — 1) only. The action is useful in our derivation 
of the perfect equations of motion. Yet we want to stress that the existence of an 
action is not really necessary. This is crucial because there is no action principle for 
the full Navier-Stokes problem. In fact, writing the above action was possible only 
because we have removed the time dependent dissipative element of the full problem 
and have forced the flow to be static. 
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What is the perfect discretization of the Poisson equation? To answer that 
question we consider the large scale lattice variable 



Vd,. = ^ / d'-'y Vdiv). (3.4) 

fd,x 

Note that the face fd^x is nothing more than a (rf — l)-dimensional cube, and Vd^x 
behaves practically as a (ci — l)-dimensional scalar. We want to impose eq. (|3.4| ) as 
a constraint for each face fd,x- For this purpose we add it to the action using a 
Lagrange multiplier field rjd^x that lives on the face fd^x 

f 1 Sp 

S[vd,Vd,Vd] = d'^'^y[-diVd{y)diVd{y) + —Vdiy)] 

J A Z Up 

+ J:^'''[vMx-^J d'-'yvd{y))-^vU (3-5) 

We have also included a term quadratic in the Lagrange multiplier field. This would 
not really be necessary, but its coefficient a will allow us to optimize the perfect 
equation's locality, and therefore it is of great practical importance. Next we derive 
equations of motion from the variational principle of the action. Varying the action 
with respect to Vd{y) yields 

didiVdiy) \-y2vd,x0d,x{y) = 0' (3-6) 

I'P X 

where Od,x{y) = 1 for ?/ G fd,x and zero otherwise. Varying the action with respect 
to rjd^x one obtains 

Vd,x - ^ / d'^-^y Vd{y) - a^arid,x = 0. (3.7) 

0- fd,x 



It appears that the constraint of eq. (|3.4|) is correctly implemented only for a = 0. 
However, varying the action with respect to Vd^x gives 

Vd,x = 0. (3.8) 

To solve the above equations we now go to momentum space. Since the velocity 
vanishes at the boundary the appropriate Fourier transform takes the form 

. d-l 

Vd{k) = / d^-^y Vd{y) n ^Hhyi). (3.9) 

Here the 

ki = —rrii, (3.10) 
Ij 

are discrete momenta with G N. Note that Vd{k) extends naturally to negative 
momentum values because 

Vd{ki, -ki, kd-i) = -Vd{ki, ki, kd-i). (3.11) 
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Similarly, for the lattice velocity field we have 

d-l 

Vdik) = a"-i Vd,. n Mk^y^). (3.12) 

X 2=1 

Now the sum is finite and extends over = 1, 2, only, where = L/a. The 
Lagrange multiplier field r]d,x is transformed analogously. In Fourier space the action 
takes the form 

Y d-l 



^^iT.i^'^dik) - va{k + 2T:l/a)Il,{k + 27il/a)\{{-iy^)r,d{k) 

k 

Vd{kf 



(3.13) 

The second sum over k is restricted to rrii G {1,2, ...,N}. In the above expression 
we have introduced 

S,{k) = n (3.14) 

i=l '^i 

as the Fourier transform of the constant 1. The function 

n,ik) Jfi'-^^^ (3.15) 

i=l ^ 

results from the Fourier transform of O^^xiy)- In Fourier space eq.( |3.(i| ) takes the 
form 

- k^ va{k) - ^6L{k) + U,{k)vd{k) = 0. (3.16) 
up 

In eq.( |3.16| ) the momentum k has the form of eq.( p.lO[ ) with uii E N being an 
arbitrary integer. On the other hand, as it stands rid{k) is defined only for momenta 
with iTLi < N . For other momenta ?7d(/c) naturally extends to 

ridiki, -fcj, kd-i) = -r]d{ki, k, kd-i), 

r]d{k + 2T:ki/a) = (-l)Sd(A;). (3.17) 

Solving eq.( p.l6| ) yields 

Vd{k) = l[n,(fc)ry,(A:) - ^-^5^]. (3.18) 

up 

We plug this back into the action and obtain 

S[Vd,Vd] = -^,Y.\-\^d{k)udd{ky^Vd{k) + ^-l/\L{k)r^^{k) + Vd{k)T^d{k)i (3.19) 
a" ^ 2 up 
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where 



1 ~''^2|„„2 



^lik) = hik + 27Tl/a) nrf(A; + 27r//a)(-l)^.K(3.20) 



Minimizing with respect to the auxihary field T]a{k) yields 

rjd{k) = Udd{k){Vd{k) + ^AL{k)). (3.21) 

up 

Reinserting this into the action finally gives 

= ^Eii^rf(^)^'i'i(^)^^(^) ^^^(^)^rf(^)]- (3-22) 

a" ^ ^ Z up 

This is a perfect action for the Poisson equation in a finite volume. The resulting 
perfect equation of motion takes the form 

ujM{k)Vd{k) + ujM{k)^—AL{k) = 0. (3.23) 
up 

By construction it is clear that this equation perfectly represents continuum physics. 
It was helpful to have an action whose variation produced the perfect equation of 
motion, but the action was not really necessary. As we will see later, all we really 
need are the equations of motion. This is important for the full Navier-Stokes 
problem, for which an action does not exist. 

For practical applications eq.( p.23D needs to be transformed back to coordinate 
space. It should be noted that ui^dik) contains the arbitrary parameter a. Hence, 
eq.( |3.23|) represents a whole family of perfect equations of motion. Now we want to 
optimize a such that the equations become as local as possible. For this purpose we 
consider d = 2. Then the sums in eq. ( |3.2CI|) can be performed analytically, and one 
obtains 

1 1 

^^22(^1)""^ = a^[ ^ 2/ /ON ~ g + ^] 



'4:Sm\aki/2) 6 

^^(^1) = "" X ■ • A /oJ (l-(-ir), (3.24) 

8sm (a/ci/2) 6 2sm(aA;i/2) 

where again ki = irmi/L. 

For the choice of a = 1/6 this gives a;22(/ci) = sin^(a/ci/2), which corre- 
sponds to the standard ultralocal second derivative in coordinate space. Transform- 
ing the whole perfect equation of motion eq.( ^.23 ) to coordinate space one obtains 

\[V2,.,+a - 2V2,^, + V2,„] = — , (3.25) 

up 
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^(1/2,1/2) 


X = 1/2 


X = 3/2 


X = 5/2 


X = 7/2 


X = 9/2 


y=l/2 


4.28599 


-0.42637 


-0.00118 


-0.00013 


-0.00000 


y = ?,/2 


-0.42637 


-0.18734 


-0.00227 


-0.00003 


-0.00000 


y = b/2 


-0.00118 


-0.00227 


0.00151 


0.00006 


0.00000 


y = l/2 


-0.00013 


-0.00003 


0.00006 


-0.00001 


0.00000 


l/ = 9/2 


-0.00000 


-0.00000 


0.00000 


0.00000 


0.00000 



Table 1: Coupling of a corner field 



for xi = {n — l/2)a with n ^ 1,N, i.e. for points away from the boundary. At the 
left boundary one finds 

^[^2,3a/2-3F2,a/2] = ^^, (3.26) 

and there is an analogous expression at the right boundary. It is straightforward to 
verify that the averaged continuum solution does in fact satisfy these equations. In 
the continuum the Hagen-Poiseuille flow in ci = 2 is given by 

V2iyi) = l%i{yi-L). (3.27) 

Consequently, the corresponding large scale lattice variable takes the form 

V2,., = - dy, V2{y,) = -A^^i - L) + -], (3.28) 

a Jxi-a/2 2 up 12 

which indeed solves eq. (|3.25|) . 

In c? = 3 the perfect equations of motion are no longer ultralocal. Still, when we 
use the optimized a = 1/6 from d = 2 also in 3 dimensions the equations turn out 
to be extremely local. In coordinate space the equation of motion takes the form 

J2ujMix,x')Vd^,, = -— Ai(x), (3.29) 

where u!dd{x,x') and Al{x) are Fourier transforms of the corresponding quantities 
in momentum space. Their values for a 3-d system with = 5 , a = 1 are given in 
tables m, 1^ and ^ for various values of x and plots for the same values are shown in 
figures 2, |^, and 4. Also the values for Dl{x) (with dp/vp set to 1) are provided in 
table 1^, and displayed in Fig. 5. 



4 Invariance of the lattice equations of motion un- 
der a factor 2 RGT 

We would like to show explicitly that the lattice equations of motions remain un- 
changed under a factor 2 renormalization group transformation. To do this we must 
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1^(5/2,5/2) 


X = 1/2 


X = 3/2 


X = 5/2 


X = 7/2 


X = 9/2 


y = l/2 


0.00151 


-0.00070 


-0.00188 


-0.00070 


0.00151 


2/ = 3/2 


-00070 


-0.19033 


-0.61801 


-0.19033 


-0.00070 


y = h/2 


-0.00188 


-0.61801 


3.24027 


-0.61801 


-0.00188 


y = 7/2 


-0.00070 


-0.19033 


-0.61801 


-0.19033 


-0.00070 


y = 9/2 


0.00151 


-0.00070 


-0.00188 


-0.00070 


0.00151 



Table 2: Coupling of a center field 



V(5/2S/2) 


X = 1/2 


X = 3/2 


X = 5/2 


X = 7/2 


X = 9/2 


y = 1/2 


-0.00118 


-0.42769 


3.85830 


-0.42769 


-0.00118 


y = 3/2 


-0.00227 


-0.18965 


-0.61602 


-0.18965 


-0.00227 


2/ = 5/2 


-0.00151 


-0.00070 


-0.00188 


-0.00070 


-0.00151 


y = 7/2 


0.00006 


0.00002 


-0.00011 


0.00002 


0.00006 


y = 9/2 


0.00000 


0.00000 


-0.00000 


0.00000 


0.00000 



Table 3: Coupling of an edge field 



* 


,T = 1 /2 


X = 3/2 


X = 5/2 


X = 7/2 


X = 9/2 


y=l/2 


0.45985 


0.66567 


0.66624 


0.66567 


0.45985 


2/ = 3/2 


0.66567 


0.99806 


0.99902 


0.99806 


0.66567 


y = 5/2 


0.66624 


0.99902 


0.99894 


0.99902 


0.66624 


y = 7/2 


0.66567 


0.99806 


0.99902 


0.99806 


0.66567 


y = 9/2 


0.45985 


0.66567 


0.66624 


0.66567 


0.45985 



Table 4: The values of Dl{x) for 5p/up set to 1 
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Figure 2: The coupling of a comer field to the lattice 

first find the recursion relation, relating the equations on the coarse lattice to the 
ones on the fine lattice. To find the recursion relation, we add to the action in the 
fine lattice variables Vd^^ (with spacing a) the constraint involving the coarse lattice 
variables V^^^, (with spacing 2a), 

K^' = ^ E Va,., (4.1) 

using a Lagrange multiplier f]d,x'- Varying the action with respect to the various 
fields yields equations analogous to the ones relating the continuum fields to the 
lattice fields 

E K^ddiV: X) + —Dl{x) - ^d,x'Vd,x'{x) = 0, 

y X 
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Figure 3: The coupling of a center field to the lattice 

Vd,x' = 0. (4.2) 

Here 9d,x'{x) — 1 ior x & x' and zero otherwise and a2 is a constant to be adjusted 
to fit with our ear her choice of optimization parameter a. Fourier transforming the 
above equations gives 

Vd{k)u;dd{k) + ^^uJddik)ALik) - n(A;)%(A;) = 0, 

d-l 

Vlik') - ^d{k' + 7rm/a)Il{k' + Tim/ a) [] (-l)i - «2r/(A;') = 0, 
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Figure 4: The coupling of an edge field to the lattice 



V(k') = 0, (4.3) 
where k'^ = Tr/Lui with Ui e (1, N/2} only, and 

m=li^os{ka/2) = ^ (4.4) 

is the Fourier transform of 0'^{x). We note that in the above equations r]d{k') natu- 
rally extends to momenta ni > N/2 by 

Uk' + nrriii/a) = {-ir^rid{k'). (4.5) 
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0.5 




Figure 5: Values of Dl{x) 
Solving the first equation for Vd{k) we obtain 

Va{k) = ^u;(A:)-^fl(fc) - ^A,(fc). (4.6) 

Inserting this relation into the other two equations and comparing with eq.(|;23D 
yields recursion formulas for uJdd{k) and Ai(A;): 

uj'ddik'Y^ = ^ Uddik' + nm/ay^Il{k' + nm/a)'^ + a2 
me{o,i}'*-i 

A'^lk') = Ai(A:' + 7^m/a)^(A;' + 7^m/a)(-l)S,'"^ (4.7) 

m6{0,l}''-i 
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We now plug in the explicit form of uJdd{k) from eq. (|3.2CI| ) into the first recursion 
relation 

^'dAk')-' = 
+ 



where we have let l[ = 2li+mi and a' = 2a in the last step. Now, letting a2 = Sa^a we 
recover the familiar form for uj'^^ik'). Thus, we see that under the factor 2 recursion 
relation, oJddik) for lattice spacing a transforms into uj'^^{k') for lattice spacing 2a. It 
can be shown that /S.L{k) behaves similarly. The fact that the equations transform 
properly indeed confirms the they are perfect lattice equations of motion. 

The recursion relation also serves as a powerful tool for numerical methods be- 
cause of the locality of the couplings. For example, one cannot solve analytically the 
lattice equations for a rectangular cross section from which a smaller rectangular 
slice has been removed. However, eq.( [4.2| ) is not specific to a given channel geometry 
and hence it still applies. To find the coupings for inward corner points one need 
only choose a very fine lattice. On it specify the couplings already known, make 
ultralocal guesses for the unknown sites, and perform several recursions. Since the 
couplings are so local, if the lattice is fine enough they serve as a very good approxi- 
mation to the actual values that define the equations of motion. Thus, the averaging 
process of the recursive step will tend to converge to fixed values, thereby providing 
the necessary information for the unknown sites. In fact, we could have derived all 
the couplings in the above tables in this fashion by making ultralocal naive guesses 
for the Poissonian and then performing recursions. 

In addition, we have numerically calculated the averaged velocity field for = 6, 
a = 1 and for = 3, a = 2. This was done by first numerically transforming uj^dik) 
and Udd{k)AL{k) back into coordinate space. The values obtained were then used 
to solve eq.( |3.29| ) for the lattice velocity field. The solution was found iteratively by 
making an initial guess for the field values and then using the equation of motion 
to find the field at x in terms of the fields at x' ^ x. With each iterative step 
the velocity field at x was replaced (if different) by the value of the field as found 
through the equation of motion. The iteration was stopped once the field converged 
to stationary values. For the above cases the field values are provided in tables |^, 
^ and plotted in figures 6 and 7. We have explicitly checked that a factor 2 RGT 
on the fine solution gives the coarse lattice fields values. The fact that they satisfy 
eq. (|3.4| ) confirms that the equations of motion found are perfect. 
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y y T, ---Trn(2A;' + 27r(2/ + m)/a)2 

2 U{2{k' + TTm)f 

11(A;' + vrm)^ 
^ 1 't}Asm\a'k', + 2nl'^ ^ 



* 


X = 1/2 


X = 3/2 


X = 5/2 


y=l/2 


0.81972 


1.40674 


0.81972 


?/ = 3/2 


1.40674 


2.48867 


1.40674 


y = b/2 


0.81972 


1.40674 


0.81972 



Table 5: Velocity field values for N = 3 



* 


X = 1/2 


X = 3/2 


X = 5/2 


X = 7/2 


X = 9/2 


X = 11/2 


y = l/2 


0.31513 


0.68669 


0.83510 


0.83510 


0.68669 


0.31513 


y = 3/2 


0.68669 


1.58922 


1.97725 


1.97725 


1.58992 


0.68669 


7/ = 5/2 


0.83510 


1.97725 


2.48758 


2.48758 


1.97725 


0.83510 


y = 7/2 


0.83510 


1.97725 


2.48758 


2.48758 


1.97725 


0.83510 


y = 9/2 


0.68669 


1.58922 


1.97725 


1.97725 


1.58992 


0.68669 


y = 11/2 


0.31513 


0.68669 


0.83510 


0.83510 


0.68669 


0.31513 



Table 6: Velocity field for N = 6 



5 The perfect velocity field 



Up to now we have shown that the perfect lattice equations of motion give exact 
results for the continuum fields averaged over cubes. In many practical cases one 
is most interested in these averaged quantities. Still, in some cases one would also 
like to know the values of the continuum fields themselves. In particular, when 
one works on a very coarse lattice the averaged lattice quantities may not contain 
enough information, and one would like to extract information for continuum points. 
Fortunately, this is possible with a perfect discretization. In fact, one can reconstruct 
the continuum field from the lattice data. For this purpose we insert eq. (|3.21| ) into 
eq.( |3.18|) and obtain 

v,{k) = hu,{k)u;,,{k){\v,{k) + ^Al(A;)) - ^6L{k)], (5.1) 
k^ up up 

which gives the continuum velocity field in terms of the lattice field V^. Transforming 
back to coordinate space one can reconstruct exactly the continuum velocity field. 
When one considers the continuum velocity field at lattice points one finds 

Vd{x) = Y.Pdd{x,x')Vd^^' H \d{x). (5.2) 

For a 2-d system pdd and are ultralocal and 

Vd{x) = l{V,^,+a + 6K^,x. + K^,x.-a) " (5.3) 

8 Sup 
for X 7^ a/2, (A^ — l/2)a and for the boundary point x = a/2 

Vdia/2) = l{5Vd,a/2 + Vd,3a/2) - ^ — • (5.4) 
o o up 
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Figure 6: The velocity field for N — 3 



For continuum points which are not lattice points one can also obtain the values of 
Pdd and Xd by Fourier transformation. 



6 Perfect discretization for static slow flows 



Now we turn to general static flows with not necessarily constant pressure gradients. 
In this case the Navier-Stokes equation reduces to 

djdjVi{y) = ^^dip{y). (6.1) 
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Figure 7: The velocity field for N = 6 

For reasons of simplicity we restrict ourselves to an infinite volume. The inclusion of 
boundary conditions is straightforward (though tedious) and can be done in analogy 
to the Hagen-Poiseuille flow. Although the above equation still follows from an 
action principle, here we restrict ourselves to working with the equations of motion 
alone. Eq.( |3.ti| ) then takes the form 

didiVj{y) - —djp + J2 Vj,xOj,x{y) = 0, (6.2) 

X 

and we still have 



20 



Similarly, the continuity equation leads to 



—dMy) + J2r)My)-0, (6.4) 



and we also have 



Px--J d% p{y) - a^Prj, = 0, Vx = 0. (6.5) 

Note that we have included a parameter /3 analogous to a to be optimized for locality 
later. It will turn out that this is not necessary, because the equations of motion 
can be optimized for locality without it. Since we are now in an infinite volume 

/d 
Vjiv) n cwi-^kiVi), (6.6) 
i=i 

/d 
d'^y P{y) n ex.p{-ikiyi), (6.7) 
1=1 

where A; e R*^. Hence, it is convenient to now let 

Vj{k) = a''-'J2V,,,Il^M-^hx^), (6.8) 

d 

P{k) = a'^J^PxU^M-ikiXi), (6.9) 

X i=l 

where k e [— 7r/a, vr/a]'^. We note that given the geometry of Fig. 1, if lie at 
points X whose coordinates are integers, then Vj^^ will lie at a point x whose j-th 
coordinate is a half integer. Also, as before, it is possible to extend the lattice fields 
to momenta outside [— 7r/a,7r/a] by letting 

Vj {k + 2nlii/a) = Vj{k) .i^j 
V^{k-r2^lf3la) = {-lfV^{k), 

P{k + 27rlii/a) = P(A;), for all i. (6.10) 
We are thus ready to go to Fourier space where the above equations take the form 

-k'^Vjik) - -^kjpik) + Ilj{k)rij{k) = (6.11) 
Vj{k) - ^Ak + 2-Kl/a)Ilj{k + 2T:l/a){-lf - a^arjj = (6.12) 

-^kjVj{k) + U{k)ri{k) = (6.13) 
P(k) -Yp{k + 27rl/a)U(k + 27rl/a) - a'^(3r](k) = 0. (6.14) 
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Solving for the continuum fields we obtain 



pik) = -{upfUikUk) - '-^k,U,ik)r^,{k) (6.15) 



v.ik) = '^kMk)vik) - ^^^^jP^U,{k)v,ik). (6.16) 

Reinserting these back into eq. ( |6.12| ) and eq. ( |6.14| ) gives equations of motion for the 
lattice fields and the Lagrange multiplier fields 

Vi{k) - iv puji{k)r]{k) + U0ij{k)r]j{k) = 

P{k) + {upfuj{k)r]{k) + iupujj{k)rij{k) = 0. (6.17) 
Here we have introduced the functions 

uj{k) = J2U{k + 27Tl/af + (3 = 1 + f3 (6.18) 

k- -\- 2Tr] ■ /n 

uJ.ik) = Y: ^k + 2nl/ay ^'^^ + 2nl/a)Uik + 27r//a)(-l)'' (6.19) 

xUi{k + 27il/a)Uj{k + 27r//a)(-l)'»+'^ + a^a. (6.20) 
We can rewrite the above equations in a simple fashion by letting Vq = P and t]q = rj 

V, = (6.21) 

where Qfj_i, is defined through the above relations. Finally, we solve for the Lagrange 
multiplier fields and set them equal to zero to obtain the equations of motion for 
the lattice fields 

n^uV, = 0. (6.22) 

We must now optimize for locality. As before, we go to = 1 and choose a and (3 
so that the equations of motion take an ultralocal form. In one dimension eq.( |6.22D 
becomes 

I ( « z/(2sin(A;/2)) \( P\^^ 

(1 + /?)«- 1/(2 sin(fc/2))2 1^ -V(2sin(A;/2) 1 + /? )\v) ^•^"'^'^^ 

Setting = and /3 = 0, we obtain familiar ultralocal equations 

2sm{k/2)V{k) = (6.24) 
2zsin(A;/2)P(A;) = (2 sin(A;/2)) V(A;). (6.25) 

Having solved this simplified version of the Navier-Stokes equations we proceed to 
add the nonlinear term. 
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7 Perturbative treatment of the nonlinear term 



Since the nonlinear term of the Navier-Stokes equation is quadratic in the velocity, 
it is possible to solve the lattice equations of motion only perturbatively. Here we 
demonstrate how to do this for the first order correction. We begin by introducing 
the modified version of eq.( |6.11| ), now with the new term 

- ^Vj{k)kjV,{k) - k'^Viik) - -^kiViik) + Uj{k)r]j{k) = 0. (7.1) 

Because we are only interested in the first order correction we need only keep terms 
in the nonlinear contribution that are quadratic in the Lagrange multiplier fields. 
Thus, for the nonlinear term, it is sufficient to express the continuum fields as linear 
functions of the Lagrange fields. We make the following ansatz for p(^) and f l^"* , the 
second order contributions: 

= Ar]{kf + r]{k)B,{k)r],{k) + r],{k)a,{k)r],{k) (7.2) 
vP = D,{k)r,{kf + r]{k)E,^7]^{k) + Fi,kVAk)Vk{k). (7.3) 

Then p*^^^ and vf* must obey the following equations, 

'-^hIi{k)W - /(^-^^■^-^''^-^■) n(fc)n,(fc)r/, - k\f^ - ^hp^^ = 0.(7.4) 

kivl^\k) = 0. (7.5) 
Solving them, we obtain 

k'' 

vf\k) = pr^(A:)^^t^n,(fc)n(fc)r^,(fc). (7.7) 

This procedure can now be repeated for the lattice fields. We make the ansatz, 

V?\k) = V,{k)A'^^{k)K{k), (7.8) 

where again Vq = P. Keeping only the second order terms, yields 

-iiypuj,V^A%V, + uJi,{k)V^Aj^,V, - pV^nl^Si^n.^V, = (7.9) 
{upfujV.A^V, + iupu,V,A^^,V, - iy^p^V^nl,Sno,V, = 0. (7.10) 

Here, the following functions have been introduced 



sW = E^^^^n.(* + 2.;/a)3 (7,11) 



, . _ ^ {kj + 27rli/a){kj + 27clj/a) - {k + 2Til/af5ij 

{k + 2.l/af 

xlii{k + 27r//a)Hj(fc + 2Til/a)Ii{k + 27r//a)(-l)''+'^ (7.12) 
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We rewrite the above equations in a more compact form 




(7.13) 



This, finally gives the first order nonlinear equation of motion 



r^M-K + n^uVxMl^V, = 0. 



(7.14) 



In principal this procedure could be repeated to obtain non-linear lattice equations 
up to any order in the lattice fields. 



8 Conclusions 

In this paper we have performed the first steps towards a perfect discretization of 
the Navier Stokes equations. The essential idea is to define coarse grained variables 
by averaging the continuum pressure and velocity fields over appropriate spatial re- 
gions. The corresponding averaged field variables naturally live on a lattice. In the 
cases of Hagen-Poisseille fiow and slow static fiows in general we derived the exact 
equations of motion for the lattice variables, and hence a perfect discretization of 
the corresponding continuum equations. The incorporation of boundary conditions 
is nontrivial, and has been done explicitly for Hagen-Poisseille fiow. For practical 
applications it is essential that the perfectly discretized equations of motion are 
very local (although in general not ultralocal). In particular, the lattice couplings 
to distant variables decay exponentially, and one can safely restrict oneself to a few 
nearby neighbors. In fact, we have performed numerical simulations using the per- 
fectly discretized equations of motion for laminar fiow in a channel with quadratic 
cross section and we have verified explicitly that an arbitrarily coarse lattice indeed 
gives continuum answers. Our method can be directly applied to more general ge- 
ometries, although non-cubic lattices require nontrivial modifications. An important 
aspect of our construction are the so-called perfect fields. Indeed our method allows 
to interpolate the continuum fields from the lattice data. Hence, especially when we 
are working on a coarse lattice, we are not limited to lattice points. 

Of course, a method like that would be most welcome in numerical simulations 
of more complicated fiuid dynamics systems, especially in the case of turbulence. 
In these simulations the effects of a finite lattice spacing are the main source of 
systematic errors. We have made a first step in this direction by deriving the per- 
fect equations of motion in the presence of the nonlinear term that gives rise to 
turbulence. However, we have included that term only to first order, while a perfect 
discretization of turbulent fiows would require a treatment to all orders. This could 
be achieved numerically in analogy to by solving a minimization problem on a 
multigrid. When one wants to simulate turbulence one should not restrict oneself 
to static flows. Then the definition of the averaged lattice variables should be gen- 
eralized. In particular, one should also average pressure and velocity in time. In 
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addition, one may want to simulate compressible fluids, which would bring in an 
averaged density variable. Our techniques are still applicable to these cases, but the 
actual implementation is nontrivial. 

In conclusion, we have proposed a renormalization group treatment of problems 
in fluid dynamics. The renormalization group implies that perfect discretizations of 
the continuum Navier-Stokes equations exist. The perfectly discretized equations of 
motion have been constructed explicitly for slow static flows, and have been opti- 
mized for locality. This is essential in practical applications. The renormalization 
group offers a systematic way of treating small scale effects in fluid dynamics. Re- 
alizing this program for the full Navier-Stokes problem requires a lot of work. Still, 
the flrst steps that we have taken, have led to promising results. 
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